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-1  In  any  inertial  navigation  system  the  platform  must  be  initially  aligned  in 
some  known  frame  of  reference  prior  to  operation  in  the  navigation  mode.  Self¬ 
alignment  methods  are  preferred  in  most  applications,  and  the  usual  procedure  is 
to  align  the  platform  locally  level  with  one  accelerometer  axis  pointing  north. 

In  the  current  generation  of  aircraft  inertial  systems  the  sensed  acceleration  is 
in  the  form  of  incremental  velocity  pulses.  These  are  accumulated  to  yield  a  measure 
of  total  velocity  within  the  granularity  of  the  incremental  velocity  change.  Kalman 
filtering  techniques  have  recently  been  applied  to  the  alignment  problem,  the  usual 
mode  of  operation  samples  the  total  velocity  at  a  fixed  rate,  and  the  resulting 
sequence  of  samples  becomes  the  input  to  the  Kalman  filter.  Granularity  in  the 
velocity  measurement  is  then  treated  as  uncorrelated  measurement  noise. 

'  VA  new  approach  to  the  alignment  problem  is  considered  in  this  report  whereby 
the  incremental  velocity  pulses  are  modeled  directly  as  the  measurement  sequence. 

This  leads  to  three  important  changes  in  the  filter  model:  (13  aperiodic  sampling 
is  obtained;  (2)  measurement  noise  due  to  granlarity  is  eliminated;  and  i^i)  a 
delayed  scate  appears  in  the  measurement  equation.  This  latter  condition  forces  the 
use  of  a  modified  form  of  the  Kalman  recursive  equations.  Results  of  Monte  Carlo 
simulations  tor  one  set  of  noise  parameters  are  given.  These  indicate  that  con¬ 
siderable  improvement  in  performance  may  be  expected  from  this  technique  relative 
to  the  more  conventional  mode. ing  technique. 
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In  any  inertial  navigation  system  the  platform  must  be  initially 
aligned  in  some  known  frame  of  reference  prior  to  operation  in  the 
navigation  mode.  Self-alignment  methods  are  preferred  in  most  appli¬ 
cations,  and  the  usual  procedure  is  to  align  the  platform  locally  level 
with  one  accelerometer  exis  pointing  north.  In  the  current  generation 
of  aircraft  inertial  systems  the  sensed  acceleration  is  in  the  form 
of  incremental  velocity  pulses.  These  are  accumulated  to  yield  a 
measure  of  total  velocity  within  the  granularity  of  the  incremental 
velocity  change.  Kalman  filtering  techniques  have  recently  been  applied 
to  the  alignment  problem,  the  usual  mode  of  operation  samples  the  total 
velocity  at  a  fixed  rate,  and  the  resulting  st 'uence  of  samples  becomes 
the  input  to  the  Kalman  filter.  Granularity  in  the  velocity  uea^ure- 
ment  is  then  treated  as  uncorrelated  measurement  noise. 

A  new  approach  to  the  alignment  problem  is  considered  in  this 
report  whereby  the  incremental  velocity  pulses  are  modeled  direct!/  3s 
the  measurement  sequence.  This  leads  to  three  important  changes  in 
the  filter  model:  (1)  aperiodic  sampling  is  obtained:  (2)  measurement 
noise  due  to  granularity  is  eliminated;  and  (3)  a  delayed  state  appears 
in  the  measurement  equation.  This  latter  condition  forces  the  use  of 
a  modified  form  of  the  Kalman  recursive  equations.  Results  of  Monte 
Carl j  simulations  for  one  set  of  noise  parameters  are  given.  These  indi¬ 
cate  that  considerable  improvement  in  performance  may  be  expected  from 
this  technique  relative  to  the  mere  conventional  modeling  teennique. 
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I.  INTRODUCTION 

Before  an  inertial  navigation  system  can  be  used  as  a  reliable 
navigation  aid,  certain  quantities  oust  be  initialized.  These  include 
the  vehicle  position  and  velocity  and  the  platform  orientation  with 
respect  to  the  navigation  coordinate  system.  Obviously,  the  accuracy 
tc  which  these  parameters  can  be  determined  is  a  limiting  factor  in 
system  performance.  If  the  vehicle  is  stationary,  we  presumably  know 
its  position  and  velocity  essentially  perfectly.  Thus,  these  values  can 
be  accounted  for  in  a  straightforward  fashion. 

To  align  the  platform  in  any  cooroinate  frame,  knowledge  oust  be 
obtainec.  of  three  angles,  usually  those  along  the  axes  of  the  given 
coordinate  system.  For  simplicity  an  x-axis  north,  y-axis  west,  and 
a -axis  up  frame  of  reference  is  used  throughout  this  paper.  Various 
leveling  and  gyrocompassing  schemes  are  well  documented  L 1_  and  will  not 
be  discussed  here.  Essentially,  leveling  is  performed  to  align  the 
platform  with  the  east-west  and  north-south  axes  while  gyrocompassing 
reduces  the  azimuth  or  "north-pointing”  error. 

Traditionally,  the  azimuth  misalignment  following  an  alignment 
period  has  been  about  an  order  of  magnitude  greater  than  the  corresponding 
level  error  angles.  In  recent  •-'ears,  however,  the  use  of  Kalman  filter 
theory  in  tstimating  the  above  quantities  has  improved  initial  alignment 
and  subsequently  increased  inertial  navigation  accuracy  capability. 

If  a  Kalman  filter  is  used  to  estimate  the  platform  orientation 
angles,  the  estimates  can  be  used  in  either  closed  or  open  loop  form. 

The  latter  method  will  be  assumed  and  is  illustrated  in  Figure  1. 


ion 
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Figure  1.  Open  loop  estimation 


In  this  scheme  output  from  the  inertial  navigation  instrument  cluster  is 
processed  in  the  Kalman  filter  for  some  predetermined  alignment  time.  At 
the  end  cf  the  alignment  interval,  platform  correction  information  pro¬ 
portional  to  the  misalignment  angles  is  used  in  the  torruer  motors  to 
correct  any  platform  misorientation. 

In  this  paper  a  gyrocoopassing  technique  using  a  new  and  unique 
incremental  velocity  measurement  algorithm  in  a  Kalman  filter  estimation 
scheme  is  presented  and  compared  with  a  system  employing  a  somewhat 
standard  time  interval  measurement  method  with  simple  periodic  sampling. 


II.  KALMAN  RECURSIVE  EQUATIONS 


Two  sets  of  Kalman  filter  equations  are  utilized.  The  first  includes 
tht  "standard”  recursive  relationships  [ 2”}  in  which  the  state  and  measure¬ 
ment  models  are  required  to  have  the  format  specified  below.  The  system 
dynamics  are  assumed  to  satisfy  the  first-order  differential  equation 

x  *  Ax  +  Df  (1) 

This  lirst-order  equation  is  then  assumed  to  transfer  to  a  discrete  state 
equation  of  the  form 

x  =  i  x  +  g  (2) 

n+1  n  n  n 

where  x  js  the  system  state  at  time  t  ,  4  the  state  transition  matrix 
n  n  n 

at  time  t  ,  and  g  ar.  uncorrelated  sequence  representing  the  system 
n  n 

response  to  white  noise  inputs.  As  in  (2),  matrix  and  vector  quantities 
will  be  implied  as  the  situation  dictates. 

The  measurement  equation  is  formulated  as 


y  =  M  x  +  6y 
Jn  n  n  n 


where  y  is  the  measurement  at  time  t  ,  M  the  linear  connection  -jstrix 
n  n  n 

relating  the  various  states,  and  dy^  an  uncorrelsted  sequence  of  measure¬ 
ment  noises  (white). 

Wi_h  the  system  modeled  by  (2)  and  (3)  the  recursive  Kalman  equations 


*  T  *  T  -1 

b  =  PM(M?H  +V)i 
n  n  n  n  n  n  n 


A  A  A 

X  =  x  +  fc  (y  -  y‘) 
n  n  n  ■'n  7n 
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*  *  T  T 

P  =  P  -  b  CM  P  M  +  V  )b 
n  n  nnnn  nn 


x*  -  i  x 
n+1  c  n 


P*  «  i  ?  i7  +  H 
n+1  nan  n 


where 


i>  =  transition  matrix  for  interval  from  t  to  t 
n  n  n+1 


x  =  tree  state  at  time  t 
n  n 

A 

xr  =  optimum  a  posteriori  estimate  of  x  at  time  t 

A 

x  =  optimum  a  priori  estimate  of  x  at  time  t 

c>n  =  gain  matrix 

v  =  measurement  at  time  t 
'  n  n 

y"  =  X  x 
J  n  n  n 

A 

?  =  covariance  matrix  pf  the  estimation  error  (x1  -  x  ) 

n  n  n 

A 

P  =  covariance  matrix  of  the  estimation  error  (x  -  x  ) 
n  an 

X  =  measurement  matrix 
n 

V  =  covariance  matrix  of  the  measurement  error,  i.e., 
n 

V„  -  E<5V*n> 

=  covariance  matrix  of  the  response  of  the  states  to  all 

T 

white  noise  drivirg  functions,  i.e.,  ~  E(gngn> 

L t  =  time  increment  between  t  and  t  x. 

n  n+1 

The  other  set  of  recursive  equations  used  involve  the  use  of  a 
modified  (delayed-state)  measurement  model.  The  measurement  equation  is 


y  H  x  +  N  x  ,+fcy 
n  n  n  n-1  n 


(9) 


That  is,  ths  measurement  at  time  t  consists  of  a  linear  connection  to 
*  n 

states  at  both  t  and  t  , .  As  before,  the  system  is  described  by  (2). 
n  n-1 

In  this  configuration  the  recursive  relationships  become  !_3j 


(10) 

(U) 

(12) 

(13) 

(14) 


A  derivation  of  the  above  delayed-state  equations  in  included  in  the 
Appendix  of  [4_: . 


! 

1 
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III.  MATHEMATICAL  MODELS 

A.  System  DyoacJ.cs  and  Incremental  Velocity 
Measurement  Model 

The  following  development  of  the  system  dynamics  is  essentially 
included  in  Pitman  [ 1 j  while  the  incorporation  of  the  new  incremental 
velocity  measurement  algorithm  is  taken  largely  from  unpublished  notes 
by  Dr.  R.  G.  Brcwn  of  Iowa  State  University. 

The  basic  equations  that  describe  the  inertial  system  error  propaga¬ 
tion  for  a  stationary  vehicle  are 
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R  =  radius  vector  from  the  center  of  the 


earth  to  ..rue  vehicle  position 


l =  local  gravity  vector,  g,  divided  by  the 


earth  radius  vector,  R 


The  above  angles  are  illustrated  in  Figure  0 


2*  z  ’  2 
(up) 


x,  y,  z  =  Ideal  coordinate 
system 


x’,  y  ,  z  =  Actual  platform 
coordinates 

x*,  y*>  2*  =  Computer  platform 
coordinates 


y  (west) 


(north) 


(:\ 

Xv ! 


Figure  2.  Coordinate  system  relationships 


It  should  be  noted  that  a  third  Schuler  equation  involving  68^  is  not 

included  because  the  vertical  error  is  assumed  to  be  zero. 

The  physical  quantities  that  we  desire  to  estimate  in  our  Kalman 

filter  are  the  level  tilts  (6  and  6  )  and  azimuth  error  (6  ).  These 

x  y  z 

angles  can  be  inserted  into  the  system  equations  through  the  relation 
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^  +  68.  The  effect  of  knowing  position  perfectly  is  that  69  =  0. 

Thus,  4  =  and  (17),  (18),  and  (19)  become 


(22) 


4  +  n  4  -  n  4 

y  2  X  X  2 


(23) 


4  +  »1  4  —  e 

2  x  y 


For  the  purposes  of  this  study  ve  make  several  simplifying  assump¬ 
tions.  tfe  first  assume  that  some  sort  of  coarse  alignment  has  taken 
place  such  that  the  level  tilts  are  relatively  small  as  compared  with 
a2imuth  error  (statistically,  at  least).  Also  we  1et  the  gyro  drifts  and 
accelerometer  biases  be  zero  (i.e. ,  we  assume  perfect  instruments). 

Then  (22),  (23),  and  (24)  become 


(25) 


(26) 


(27) 


An  explicit  solution  can  now  be  written  for  4  ,  4  ,  and  4  . 

x  y  z 


o 

O 

V-' 

•o* 

II 

✓—s 

■** 

(28) 

** 

N 

/-N 

rt 

N 

N 

O 

(29) 

t 

*  (t)  =  '  a  *>  dt  +  4  (0) 
y  %  X  z  y 

=  a  4  (0)t  +  4  (0)  (30) 

x  z  y 


or 


s 


6 

X 

10  0 

i 

<-N 

o 

zero 

N 

- 1 

o  i  n  t 

X 

0  0  1 

# 

iy(0) 

v°> 

+ 

driving 

function 

Ol) 


These,  then,  are  the  basic  state  equations. 

The  measurement  model  must  take  into  account  the  "real-world” 
accelerometer  mechanization.  Basically  the  accelerometer  emits  a  pulse 
whenever  a  predetermined  increment  in  velocity  (integrated  acceleration) 
is  reached.  Hiysicallv  this  can  take  the  form  of  an  integrator  at  the 
accele* oo.eter  output.  The  integrator  is  used  in  conjunction  with  a 
threshold  detector  that  emits  a  +  LV  pulse  whenever  the  preset  threshold 
is  reached.  When  the  delta  velocity  pulse  is  emitted,  the  appropriate 
positive  or  negative  LV  increment  is  fed  back,  ideally  resetting  the 
integrator  to  zero.  An  example  of  an  integrated  acceleration  function 
and  the  corresponding  sequence  of  LV  pulses  is  shown  in  Figure  3. 


Integrated 

Acceleration 


Accelerometer 

Ojtput 


LV 


•LV 


1 

LI _  .. 

tlt2  C3  * 

* 

Tine 

Figure  3.  Integrated  acceleration  function 
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Note  the  nonunifora  time  interval  between  pulses. 

Using  the  open-loop  correction  scheme  discussed  .ibove,  we  want  to 

estimate  the  level  tilts  (6  ,  4  )  and  the  azimuth  error  (6  )  based  on 

x  y  z 

a  sequence  of  AV  measurements  that  occur  in  a  finite  interval  of  time. 

Explicitly  the  AV  measurement  is  (for  both  accelerometers) 
t 


and 


Now 


n 


a  (t)dt 
x 


n-1 


a .(t)dt 


“r.-l 


ax  ~  ~  8*y  +  (Noise 


«  =  Z*  +  (Noise) 

>  A  y 


(32) 

(33) 


Using  (31),  (32),  and  (33)  the  problem  decouples  as  follTws 


East -Vast 
Channel 


r  6  =  i  (0) 

i  X  X 

/ 

[a  =  +  (Noise) 

y  *  y 


North-South 

Channel 


'i  n  t" 

m. 

■— N 

o 

_ J 

X 

y  i 

* 

• 

I 

-J 

-° 

N 

© 

-  g^  +  (Noise) 
v  x 


Process  Xodel 

Measurement 

Process  Model 

Measurement 


At  this  point  we  will  develop  only  the  north-south  measurement  equations 
as  that  channel  furnishes  the  information  concerning  azimuth  error.  Tne 
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east«vest  derivation  proceeds  in  an  analagous  fashion  and  will  be 
omitted . 

The  actual  measurement  is 
t 


r. 


a  dt 
x 


'n-1 


Writing  this  out  explicitly  we  obtain  (from  (32)) * 


u  u  u 

f  a  (t)dt  *  j  -  g<  (t)dt  +  J  noise 


n-1 


n-1 


n-1 


Substituting  for  rfy(t). 


u  ti 

ax(t)dt  *  -  g  j  [&x<,(0)t  +  d  (O)Jdt  +  j  noise 


n-1 


n-1 


n-1 


-  %«c  (0>(t*  -  t2  ,)  -  grf  (0)(t  -  t  -) 

x  z  n  n-1  y  n  n-1 


I 


noise 


n-1 


Now  substituting 


v°>*  Vl)  ‘°A‘ 

-  VV  '‘Wo 


(34) 


and  d  (0)  =  6  ( t  )  into  (34)  we  obtain 
z  z  n 


2  2 


J  *,<4)*  -  -  48  “A".  -  ‘»-l>  -  sW/tn)  '  -  Vl> 


o-l 


+  J  noise 
tn-l 


(35) 
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Rearranging  (35) 
t 

n 


Jt  \M6t  '  8  nA(t„  -  *  %(t„  +  Cn-1^ 

n-1  tn 


-  g  i  (t  )(t  -  t  )  *■  noise 

y  n  n  n-l  * 

ta-l 


(36) 


Sow  let  t  -  t  .  =  Lt  .  Then 
n  n-1  n 


; n  *x(t>dt  -  («  c.ya,w  *  (-  8*taxy(tn) 

Cn-1  t 

,  n 

+  noise 

at 

cn-l 


(37) 


Equation  (37)  is  now  in  the  correct  measurement  foraat  except  for  the 
noise  tern. 

To  account  properly  for  the  integrated  noise  tern  in  (37),  let  us 
assume  that  the  major  source  of  noise  is  random  lateral  notion  of  the 
vehicle  plus  white  instrument  noise.  Although  lateral  notion  was  ignored 
in  the  system  process  derivation,  it  is  an  important  noise  source  due 
to  wind  buffeting,  loading,  and  general  random  notions  that  occur  to  the 
vehicle.  In  deriving  our  noise  model  we  will  let  the  vehicle  be  an 
aircraft,  thus  knowing  that  there  is  no  appreciable  net  random  notion 
(position)  when  it  is  stationary  on  the  ground.  At  the  acceleration 
level  this  is  a  special  type  of  noise;  it  is  such  that  the  double  integral 
of  the  acceleration  noise  is  bounded,  and  is  thus  a  stationary  process. 
Therefore,  let  us  assuae  that  at  the  position  level,  the  process  is 
shaped  Markov.  Tnen,  so  that  acceleration  and  velocity  have  bounded 
variance,  we  postulate  the  model  shown  in  Figure  ~  (position  that  is 
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Markov  shaped  by  a  second-order  filter).  Kote  that  position,  velocity. 


Unity 
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Koise 

f(t) 


- 
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r 
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2  2 
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Koise 


^(t) 


Figure  4.  Koise  aodel 

2 

and  acceleration  all  have  bounded  variance  with  the  paraaeters  n  ,  r, 
ii  and  Q  chosen  to  fit  the  physical  situation  at  hand. 

The  seccad  part  of  the  integrated  noise  tern  is  the  contribution 
due  to  instrument  noise.  This  noise  is  assured  to  arise  fron  the  basic 
accelerometer  itself  (e.g. ,  proof  cass  jitter  due  to  noise  fron  the 
servo  electronics)  and  thus  appears  as  integrated  acceleration  noise 
at  the  instrument  output.  We  will  assume  that  the  basic  acceleration 
noise  is  white,  thus  the  accelerometer  output  is  a  random  walk  process 

the  physical  situation, 
first-order  differential 

equation 

n2(t)  =  v(t)  (38) 

Integrating,  we  obtain 


as  shown  in  Figure  5.  is  cho3en  to  fit 
The  random  walk  process  is  described  by  the 


n^(t)  =  ^(O)  +  h(t) (Driven  Response) 


(39  > 
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or  in  difference  equation  format 


n  =  n  +  h 
2r+1  2n  n 


where  h  is  found  by  forming  the  convolution  integral  [5] 


h  -  I  y(u)v(£t-u)du 
°  0 


y(u)  is  the  weighting  function  of  the  filter  (i.e.,  the  inverse  Laplace 
transform  of  1/s). 


Instrument  Noise 


Acceleration  Noise 


(white  noise -amplitude  K^) 


Figure  5.  Instrument  noise  random  process 


y(u)  lCl/s] 


l*v(At“u)du  =  0 


since  the  white  noise  input  v(t)  has  zero  mean.  Also, 

At  it 
2  ^  m 

=  J  j  y(u)y(w)v(it-w)v(it-u)dudw 

0  0 


it  it 


h2  =  I  y(u)y(w)v(it-u)v(it-v)dudw 

n  •  * 


Since  v(t)  is  white  noise  with  amplitude  v(t-u  )v(t-w)  =  Kji(u-w), 

where  o(  )  is  the  Dirac  delta  function,  Then 


it  At 


1* 1*K^6  (u-w)dudv; 


0  0 


Using  sifting  integral  theory. 


-  J  K,*> 


-  V' 


Thus  h^  can  he  represented  fay  an  uncorrelated  sequence  of  Gaussian  random 
variables  with  zero  mean  and  variance  equal  to  K^t. 

Let  us  now  assign  states  to  the  physical  quantities  to  obtain  the 
correct  process  model  format  for  r'.e  Ka  mar  filter.  Let 


Figure  can  be  redrawn  as 


I  M  •  a2 

i  r 

f(t)  *  2  2~ 

0  nity  ”  r  r 

White  _ _ 

Noi sc ' 


Acceleration 

Noise 


Velocity- 

Noise 


1  | — >  Position 


Fijture  6.  State  formulation  of  noise  mode 


=  z(t) 

(position  noise) 

(50) 

x4  =  Z(t) 

(velocity  noise) 

Oi) 

X 

II 

n: 

rt 

(acceleration  noise) 

(52) 

differential  equation  is 

z  +  (3  +  2£u>r)z 

2  *  2 
+  (u;  +  23£uii  )z  +  3u;  2  = 
r  r  r 

(53) 

Then  in  matrix  form  we  have 


Equation  (54)  is  in  the  standard  state  variable  form  of  Equation  (1). 

The  state  transition  matrix  can  now  be  found  in  closed  form  from 

*(t)  »  -  A)"1]  (55) 

where  I  is  the  identity  matrix.  Each  term  of  ^(t)  is  of  the  fora 

K.j  e’6At  +  L±i  e"CU;rAt  sin  Vl^C^-At  +  0^) 


(56) 


17 


where-  K.  . ,  L. . ,  and  9..  are  functions  of  a1  ,  '  r,  and  . 

XJ  lj  r*-** 

The  driven  responses  for  the  states  (g  )  can  be  found  in  closed 

v  n 

form  by  again  forming  the  convolution  integral 
t 

x^t)  =  J  y^(u)f (t-u)du  (37) 

0 

where  y.(u)  is  the  weighting  function,  between  f(t)  and  the  state  of 
interest,  x..  Then 

l 

t 

x.(t)  =  y. (u)f (t-u)du  =  0  (58) 

i  v  i 

0 

and 

_  t  t 

x^(t.>  =  r  ^  yi(u)yj/v)  f  (t-u)f (t-v)dudv  (59) 

0  0 

~2 

as  before.  x.(t)  again  indicates  taking  the  expect- d  value.  Since  f(t) 
is  unity  variance  white  noise ,  f (t-u)f (t-v)  is  just  6(u-v).  Equation  (59) 
then  be  cones 

_  t  t 

x?(t)  =  •  y.(u)y  (v)5 (u-v)dudv  (60; 

X  •*  c  1  1 

o  o 

t 

=  y?(v)dv  (61) 

v  1 

0 

Thus  the  driven  response  for  each  state,  x^,  is  a  random  variable  with 

t 

,  2 

zero  mean  and  variance  equal  to  ^  y^(v)dv,  where  we  let  t  =  it. 

0 


Recapping,  then. 


rfn(3,3) 

through 

Vs -5> 

(Transition 

Matrix) 


X4  +  *4 

a  r 


The  last  state  to  be  considered  is  described  by  (40).  Here  we  let 


n2  =  x6  and  hn  *  86  •  Th  ■ 


n  n 


*  xa  +  8, 


Finally,  combining  all  six  states  in  matrix  form. 


I  cix&t  0  0  0  0 

0  1  0  0  0  0 


0  0  |  <5q(3,3)  ,0 

.  i 

0  0  J  through  t  0 

0  0  l  i  (5,5)i  0 

1 _ n_  _  . 

0  0  0  0  0  0 


x  + 


Tne  process  model  is  now  in  the  correct  format. 

Returning  to  the  measurement  equation,  (37),  we  now  write 


y(tn) 


•*dt  "  <-  ^VXl<tn>  +  ^■Vtn)x2(tn) 


+ 


(Accelerometer  Noise )dt 


But  the  integral  of  accelerometer  noise  is  just  x^  +  x^,  thus* 

y<*;>  -  +  tV'ta),<2<tn) 


+  X4^n^  ~  x4^n-l^  +  X6^n^  ‘  x6*fcn-l^ 


Note  the  delayed  states  in  the  measurement  equation  and  the  absence  of 

the  "standard"  measurement  white  noise.  It  should  be  noted  from  (62), 

however,  that  x,(t  )  -  x,(t  , )  is  just  g, (t  ,)  or  an  uncorrelated 
on  o  n-l  o  n*l 

sequence  of  Gaussian  randos  numbers  *rith  zero  mean  and  variance  given 
by  (47).  Thus  this  quantity  could  be  accounted  for  in  the  measurement 
equation  by  denoting  it  as  Sy^  (measuretjent  noise).  Then  (64)  becomes 


y(t„)  =  <'«6tn>xi(tn)  +  (^8iVt«>x7<t„) 


x  a  z  n 


+  x4(tn>  (65) 

Including  state  six  as  a  noise  contribution  in  the  measurement  equation 
contains  a  significant  advantage  when  implementing  a  recursive  routine 
as  this  technique  reduces  the  dimensionality  of  the  problem  by  reducing 
the  order  of  the  state  process.  However,  f64)  will  be  considered  the 
measurement  equation  for  purposes  of  this  paper  because,  os  will  be 
shown  later,  the  incremental  time  measurement  noise  cannot  be  handled 


as  in  (65)  above. 

The  model  is  now  complete  for  us^  in  the  delayed-state  Kalman 
recursive  equations,  i.e.. 
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y  =  1-gAt  %g  Ci  0  1  0  1]  x 

n  n  x  n  n 

-*-1  0  0  0  -1  0  -1  }  xn  (66) 

n-i 

Also  note  that  each  measurement  is  either  a  positive  or  negative  delta 
velocity  quantum  and  that  the  measurement  time  interval,  AtQ,  is  non- 
uniform  and  varies  as  the  process  evolves. 

B.  Incremental  Time  Measurement:  Model 


As  seen  in  the  previous  section  the  incremental  velocity  algorithm 
involves  delayed  states  and  nonstandard  time  intervals.  Schemes  used 
for  Kalman  filter  gyrocompassing  in  previous  studies  have  employed  other 
measurement  models,  usually  involving  a  Aniform  time  interval.  One 
such  model  is  described  in  some  detail  in  [ 6j . 

In  this  model  the  output  pulses  frcm  the  accelerometer  are  stored 
in  a  digital  pulse  count  register.  The  register  is  then  sampled  at  some 
uniform  rate  to  furnish  the  Kalman  filter  measurement.  The  pulse  count 
in  the  register,  then,  is  a  measure  of  the  total  velocity  (integrated 


acceleration)  and  not  the  incremental  change.  The  measurement  equation  is 


y(tR)  *  dt  +  (Integrated  Acceleration  Noise)  +  n^  (67) 


The  first  term  is  handled  as  in  the  previous  section  except  that 
the  time  interval  between  measurements,  At,  is  replaced  by  the  total 
elapsed  time  t.  The  second  quantity  in  (67)  is  just  x^(tR)  +  x^(£  ) 
as  derived  in  Section  A.  Notice  that  the  instrument  noise,  x^,  must  be 


handled  as  a  state  since  It  is  a  random  walk  process  and  not  a  sequence 
of  uncorrelated  random  numbers.  The  n^  term  represents  the  quantitation 
error  at  each  sampling  of  the  accelerometer  pulse  count  register. 

Figure  7  illustrates  the  measurement  procedure  for  a  typical  member 
function. 

The  registers  are  sampled  at  uniform  At  intervals  and  it  can  be 
seen  that  the  measurement  may  be  in  error  by  as  touch  as  +  AV.  Although 
this  measurement  noise  is  not  normally  distributed  with  zero  oiean,  and 
is,  in  fact,  uniformly  distributed,  the  quantization  error  is  assumed 
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to  be  normally  distributed  l6j.  Thus  a  is  assumed  to  be  zero  *>an 

13 

Gaussian  white  noise  with  variance  (AV)2.  Letting  n  =  6y  we  now 

®  n 

write  (67)  as 

y(tn)*[-gt  %gnyt2  0  10  1]  x(tn)  +  oyn  (68) 

Equation  (68)  is  new  in  the  required  format  for  the  standard  recursive 
Kalman  equations.  6y^  is  the  measurement  noise  with  variance  (AV)2. 
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IV.  DIGITAL  SIMULATION 
A.  Method  of  Analysis 

Evaluation  of  the  oev  incremental  velocity  algorithm  using  any  sort 
of  analytic  technique  is,  at  best,  very  difficult-  For  this  reason, 
cce^>arison  vith  a  more  or  less  ^standard"  estimation  technique  was  chosen 
as  the  best  means  for  evaluating  system  performance.  The  approach  taken 
was  to  simulate  the  physical  random  process  and  use  both  estimation 
techniques  simultaneously.  At  the  end  of  same  predetermined  alignment 
time  the  ability  of  the  estimators  to  determine  the  actual  value  of  the 
azimuth  error  and  level  tilt  was  compared .  Ten  Monte  Carlo  simulations 
were  performed  to  establish  some  measure  of  statistical  validity.  The 
actual  computer  program  •  »ed  for  the  simulation  is  included  in  Appendix  A 
wile  a  flow  chart  illustrating  program  organization  and  the  various 
functional  blocks  is  shown  in  Figure  8. 

The  program  itself  consists  of  four  sections,  those  being  the  main 
program  and  three  subroutines.  The  main  progtam  simulates  the  physical 
syttea  and  the  delta  velocity  and  time  interval  measurements  while  the 
subroutines  handle  the  Kalman  estimation  and  the  computation  of  the  state 
transition  and  H  matrices  in  closed  form  (the  H  matrix  is  the  covsriance 
matrix  of  the  state  response  due  to  the  white  noise  inputs). 

B.  Process  and  Measurement  Simulation 

The  system  is  modeled  by  the  difference  equation  described  by  (2) 
and  uses  a  delta  time  interval  of  one  millisecond  (Table  1).  It  should 


24 


Figure  8.  Computer  flowchart 


be  noted  from  Appendix  A  that  x,  is  not  computed  in  a  standard  iterative 

iri-1 

Banner  f roa 

*1  -  *1  +  ‘V*t-,2  <6,) 

o4-l  n  n 

The  second  ten  in  (69)  is  very  saall  (typically  2  x  10  )  and  when  x^ 

becoaes  large  the  technique  is  inaccurate  due  to  round-off  errors.  This 
problem  could  be  handled  by  using  double  precision  techniques  or  by 
computing  the  second  term  at  each  iteration  and  adding  it  to  the  fixed 
initial  x^(0).  The  latter  approach  is  utilized  in  Appendix  A. 

For  purposes  of  the  simulation  we  assuae  that  the  accelerometer 
instrument  noise,  x^,  is  zero,  i.e. ,  again  we  assuae  perfect  instruments . 
Thus  on*y  five  states  are  siailated  in  the  program  The  state  transition 
matrix  is  computed  from  H (si  -  A)  *]  as  explained  in  Section  111. 

Noise  parameters  chosen  to  fit  the  physical  situation  are  included  in 
Table  1. 

Table  1.  Noise  parameters 


Parameter 

Value 

2 

o 

9.0  inches^ 

a 

0.1  sec  * 

c 

0.5 

IP 

r 

0.1  cps  *  0.21T  rad  lacs /second 

At 

0.001  second 

latitude 

40°  N 
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The  selected  parameters  represent  a  vehicle  placed  in  a  reasonably 
benign  environment.  A  a  example  of  this  situation  night  be  a  large 
flexible  aircraft  sitting  on  a  ran?  with  moderate  vind,  loading,  etc. 
causing  sane  randoea  motion.  The  At  increaent  is  chosen  so  that  the 
process  looses  essentially  continuous  with  respect  to  the  system  dynamics. 
Closed  form  equations  for  the  state  transition  matrix  using  the  above 
parameters  are  given  in  the  main  program  and  HPHi  subroutine  and  are 
included  in  Appendix  A. 

Since  ve  performed  Monte  Carlo  simulation,  random  initial  conditions 

vere  supplied  tc,  the  process  for  each  run.  Uncorrelated  unity  variance 

random  numbers  with  Gaussian  distribution  vere  available  on  magnetic  tape 

and  vere  transformed  via  the  Schsridc  orthogonalization  procedure  7j  into 

initial  conditions  with  the  desired  variances  and  covariances.  In  order 

to  accomplish  this,  let  z^,  z2  ,  z£»  and  z5  De  uncorrelated,  unity 

variance,  zero  rean  Gaussian  randoo  variables.  States  x^  and  x2  are 

decoupled  from  the  other  states  and  are  assumed  uncorrelated.  Their 

2  l 

assumed  variances  are  (1  tain)  and  (60  min)  respectively,  thus  ve  let 
xx(0)  *  z, 

and 

x2(0)«60z2  (70) 

States  x,,  x, ,  and  xr  are  not  uncorrelated.  Tnus  these  initial  conditions 
5  u  j 

must  reflect  the  appropriate  cross -correlation.  Let 


|  rr-1 
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x3(0) 


x4(°) 


x-(0) 


a2  a3 


bl  °2  b3 
C1  C2  C3 


Z4  I  (7l> 


with  a^,  a^,  and  set  equal  to  zero.  Then 


x3(0)  =  a^ 


2 .  .  2  2  2 
x3(0)  =  a  z,  =  a 


This,  then,  specifies  a^.  In  a  similar  manner 


x3(0)x4(0)  =  aibxz3  +  ai^*3z/»  s  aibi  (73) 

as  z3  and  z^  are  uncorrelated  and  z3  has  unity  variance.  Since  a^  is 
specified  by  (72),  can  now  be  found.  The  remaining  members  of  C  may 
be  obtained  in  a  similar  fashion. 
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Although  the  terms  in  the  g^  column  vector  were  derived  in  Section  III, 
the  use  of  some  small  time  interval  approximations  simplify  the  computa¬ 
tional  effort  (the  derivations  in  Section  III  will  be  useful  later  in 
computing  the  H  matrix).  If  the  system  is  described  by  (1),  then. 


integrating  x  yields 
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At  At  At 

^  x  dt  *  w  A(t)x(t)dt  +  ^  Df(t)dt 
0  0  0 

For  small  At, 

At 

A(t)x(t)dt  «  A(t)x(0)At 

v 

0 

In  addition,  D  is  independent  of  t,  so  (75)  becomes 

At 

x(At)  -  x(0)  as  A(t)x(0)At  +  D  '  f(t)dt 

0 


x(At)  ss  (I  +  AAt)x(O)  +  Du  (At) 


where 


At 

u(At)  *  '  f (t)dt 

w 

0 


In  difference  equation  format 


x  ss  (I  +  AAt)x  +  Du 
n+I  n  n 

To  simulate  u  we  take 
n 

_  At  At  _ 

u2(At)  =  **  '  f  (u)f  (v)dudv 

w 

0  0 
or 

_  At  At 

u*"(A»:)  *  t(u-v)dudv 

0  0 


(75) 


(76) 


(77) 


(78) 


(79) 


(80) 


(81) 


»*t+np*?ty 
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since  £(t)  is  unity  white  noise.  Evaluating  the  integral  as  before 


u2  (At)  =  ht 


(82) 


Thus  is  simulated  as  an  uncorrelated  sequence  of  normally  distributed 
randan  numbers  with  At  variance.  Comparing  (79)  with  the  difference 
equation  format,  (2),  we  see  that 


g 


n 


*  Du 

n 


where 


(83) 


(84) 


i 

i 

1 

I 

i 

! 

t 


D  appears  here  as  a  single  column  vector  with  five  entries  since,  as 
previously  stated,  we  ignore  the  accelerometer  instrument  noise. 

As  noted  in  the  flowchart  there  is  a  difference  between  the  true 
value  of  vehicle  velocity  and  the  accelerometer  register  output.  The 
accelerometer  only  starts  integrating  acceleration  at  the  instant  the 
instrument  is  engaged,  i.e.,  at  t=0.  No  knowledge  of  acceleration  prior 
to  t=C  is  available,  thus  the  accelerometer  does  not  have  a  measure  of  the 
true  velocity  initially.  The  "residual  initial  condition"  term  appears 
because  it  is  realistic  to  assume  that  the  accelerometer  integrator  is 
not  zeroed  perfectly  at  time  t*0,  i.e.,  some  initial  random  voltage 
within  the  AV  granularity  is  present  on  the  integrator  output.  The 
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relation  between  these  quantities  for  a  typical  velocity  function  is 
illustrated  in  Figure  9.  The  6V  quantization  used  for  the  simulation 
is  taken  to  be  the  sane  as  in  [6j  and  is  listed  in  Table  2. 

C.  Delta  Velocity  Kalman  Estimator 


The  subroutine  DVKAL  is  used  to  implement  the  incremental  velocity 
Kalman  estimator.  The  recursive  equations  used  are  slightly  modified 
from  those  given  in  (10)  through  (16).  The  relationships  are 


n-1 


(85) 


* 

P 

n 


b 

n 


*  T  T 

(PH  +<  .P  ,HA)Q  1 
n  n  n-1  n-1  n 


+  b  <y_ 

n  n 


(86) 

(87) 

(88) 


P  =  P*  -  b  QbT  (89) 

n  n  n  n 

with  the  quantities  as  previously  defined.  Modification  is  required 
because  the  "step-ahead"  feature  of  the  Kalman  filter  is  not  possible, 
as  the  time  increment  between  measurements  is  not  known  unt i 1  the  actual 
measurement  is  made.  Thus,  it  is  noc  until  time  t^  that  we  have  access 
to  the  value  t  -  t  The  a  priori  estimates  of  x  and  P  are 

obtained  at  t  by  evolving  the  optimum  estimates  at  t  through  (85) 


and  (86). 


Integrated 

Acceleration 


Figure  9.  True  velocity  and  accelerometer  output 

Certain  initial  quantities  are  needed  to  begin  the  recursive 

A 

routine.  As  seen  from  (85)  and  (86).  initial  values  of  xr-1  and  Pr_1 

as  veil  as  values  for  H  x  and  mist  be  obtained.  Since  the  means 

of  all  the  states  of  x  are  zero,  x  *  0  is  chosen  as  the  initial  best 

n  n-i 

estimate  of  the  states  at  t-0.  The  values  for  the  elements  of  the  initial 
error  covariance  matrix  are  chosen  to  be  the  steady-state  variances  and 
covariances  relating  the  various  states.  We  assume  that  the  initial 
"coarse"  alignment  period  has  reduced  the  level  tilt  (x^)  to  a  value 


32 


Table  2.  Quantities  used  in  the  simulation 

DVKAL 


"1  0 

0 

0  0 

0  3600 

0 

0  0 

(Initial)  = 

0  0 

8.808 

0  -.47741 

0  0 

0 

.47741  0 

0  0 

-.47741 

0  .2184 

Pn-l(1*l) 

in  (uinutes 

of  arc)"" 

Pn-1(2’2> 

in  (minutes 

of  arc)2 

P  ,(3,3)  in  (inches)2 
n-1 

2 

P  ,(4,4)  in  (inches/second) 
n-1 

?  ,(5,5)  in  (inches /second2) 2 

n-1 

V  =  E(cy  5yT)  =  (.1AV)2  =  (.0394)2 
n  n  n 

M  =  (-git  %g  Q  At2  0  1  0  ) 

n  n  x  n 

N  =  (0  0  0  -1  0  ) 

n 

A 

x  =  0 
n-1 

£V  =  1  cm/sec  =  (1/2.54)  in/sec 


_ TIKAL 

* 

P  (Initial)  =  Sane  as  P  ,  above 
n  n-1 

V  *  (£V)2  *  (1/2.54)2  in2 /sec2 
n 

M  =  (-gt  -hg  t2  0  1  0  ) 

n  x 

At  =  1  second 


33 


having  a  standard  deviation  of  one  minute  of  arc  while  the  aziaith 

error  (x£)  is  assumed  to  have  a  standard  deviation  of  sixty  minutes  of 

arc.  Although  we  assume  zero  initial  cross -correlation  between  x^  and 

x^,  we  oust  take  cross-correlations  of  the  noise  states  into  account  as 

the  white  noise  driving  function  distributes  itself  intu  these  states 

as  the  process  evolves.  Values  for  these  steady-state  noise  variances 

and  covariances  may  be  obtained  by  choosing  some  arbitrary  ?n  .  matrix 

and  letting  (86)  evolve  until  the  values  converge.  Note,  however,  that 

as  time  approaches  infinity  approaches  zero,  thus  PR(t  —  00 )  is 

equal  to  H  (t  —  <=).  These  values  are  readily  confuted  via  the  HHil 
n 

subroutine  (discussed  in  E  below)  and  are  given  in  Table  2.  The  H  . 

n-i 

and  4  .  matrices  are  obtained  from  the  HEHI  subroutine  each  time  the 

n-1 

equations  are  used.  Other  quantities  needed  for  the  estimation  are 

V  ,  M  ,  and  N  and  are  included  in  Table  2.  As  noted  in  the  table, 

n  n  n 

V  is  not  tero  as  was  assumed  in  Section  III.  The  desire  to  insure  a 
n 

"safe"  measurement  model  for  the  simulation  led  to  the  addition  of  a 

small  amount  of  measurement  noise.  If  the  measurement  noise  is  zero 

we  see  from  (15)  that  as  P  ,  decreases,  Q  becomes  very  small.  Since  x. 

n-1  1 

and  x^  are  deterministic,  these  elements  of  eventually  approach 

zero  making  the  value  of  Q  subject  to  computer  round-off  error.  Q  can 
in  fact  become  negative,  causing  the  Kalman  estimation  to  diverge.  In 
addition  to  the  divergence  problem,  almost  any  physical  situation  has 
some  small  amount  of  random  noise  nresent  at  the  measurement.  The 
addition  of  the  measurement  noise  thus  makes  the  measurement  model 
"safer"  in  a  physical  sense.  The  magnitude  of  the  noise  is  assumpd 
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to  >e  10%  of  the  inciernental  granularity. 

Since  we  assume  that  the  accelerometer  integrator  ia  not  perfectly 
zeroed  at  t«0,  the  first  AV  meaourement  must  be  ignored  as  the  At  time 
interval  at  the  time  of  the  measurement  does  not  correspond  to  a  true 
velocity  increment. 

D.  Time  Interval  Kalman  Estimator 

The  TIKAL  subroutine  is  included  in  Appendix  A  and  implements  (4) 

through  (81  to  estimate  the  various  states.  The  accelerometer  pulse 

count  register  is  assumed  to  be  sampled  every  second  [t>],  (i.e.,  At  =  1 

second),  thus  H  and  i  are  constant  matrices  with  values  obtained  from 
n  n 

the  HPHI  subroutine.  The  initial  a  priori  estimate  of  x  is  chosen  to 

n 

★ 

be  zero  while  the  initial  error  covariance  matrix,  P  ,  is  obtained  as 
in  Part  C  above.  These  and  other  quantities  used  in  the  estimation  are 
given  in  Table  2. 

E.  HPHI  Subroutine 

The  Hihl  subroutine  computes  the  state  transition  (PHI)  and  driven 
response  covariance  (H)  matrices  in  closed  form,  given  a  specified  time 
increment.  As  previously  stated,  KJI  is  obtained  from  the  inverse 
Laplace  transform  of  (si  -  A)  *  and  in  general  form  is 


1 

n  At 

X 

0 

0  0 

0 

1 

0 

0  0 

0 

0 

f  A 

1  *33 

(At) 

0 

0 

1 

through 

0 

L 

0 

J 

» 

*55  (At) 

*(At)  = 


(90) 


nerK’*ir~TriJ^mJ‘J-LZ’JL£ri±iJCMK3*~. 
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As  noted  in  Section  III  i^(C^.t)  through  <S^(At)  are  of  the  form  (56) 

where  K„  ,  ,  and  8^  have  been  previously  computed  and  are  functions 

2 

of  £,  3,  u)r,  and  a  . 

T 

The  H  matrix  is  defined  as  E(g  g  ),  the  covariance  matrix  of  the 
n  °n  cn 

state  response  due  to  the  white  noise  input.  The  first  two  states  are 
completely  decoupled  from  the  white  noise  input ,  thus  the  corresponding 
elements  in  the  H  matrix  are  zero.  The  variances  and  covariances  relating 
the  final  three  states  are  obtained  by  forming  convolution  integrals  in 
the  manner  of  (57)  through  (61).  In  general. 


Ac  At 

8igj  =  XiXj  =  I  J  yi(u)yj(v)6(u-v)dudv 
0  0 


At 


“  t 

O 


yi*V^yj^dv 


(91) 


(92) 


0 


where  y^(u)  is  the  weighting  function  relating  f(t)  and  x^.  The  various 
weighting  functions  are  obtained  from  Figure  6  with  (92)  having  the 
g- neral  form 


-  ?  4' 

8igj  =  3ur 


Kij  6 


-23  At 


+  L  e 
ij 


. 


-(3+Cu;r)At  r-(3+C«ur)  sin  (DAt+e^  -  u)r  Vl~C  cos  (DAt+8^ 


+  M.,  .  e 
ij 


-(3-KuJr)At  r-(3+Car)  sin  (DAt+6^ )  -  u>r  71  -f2  cos  (DAt+8^) 

C 
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’■•fjpjii&rymJWl*  .s*  w*  -,' 


.  -2Coi  6t 

+  N£j  co*  (ei  -  ej)  e  r 


+  p.  (&**.».>  ,1 

- - ^~2 - - - ( 

1  «■£  J  I 

J  0 


where 


c  e  52  +  2£u>r&  +  «j2 


D  *  -rlC? 

The  constants  K.  .  t  m  D 

ij  ’  ij’  ij’  ij»  0i»  and  ®j  are  again  precomputed 
and  are  functions  of  C,  3,  »  ,  and  u2. 


37 


V.  RESULTS 

The  performance  of  the  incremental  velocity  Kalman  estimator  (DVKAL) 
was  evaluated  by  comparison  with  the  "standard"  time  interval  (TIKAL) 
system.  The  results  of  the  ten  Monte  Carlo  simulations  are  shown  graph¬ 
ically  in  Appendix  B.  The  DVKAL  and  TIKAL  estimates  as  well  as  the 
true  values  of  states  one  and  two,  are  plotted  versus  time  in  these 
graphs.  In  addition  to  the  state  estimates,  of  course,  the  Kalman  filter 
also  computes  the  error  covariance  matrix.  These  error  covariance  terms 
for  level  tilt  and  azimuth  error  are  averaged  over  the  ten  runs  and 
plotted  with  the  actual  squared  error  in  Figures  10  through  13. 

TWo  factors  oust  be  considered  when  evaluating  system  performance, 

these  being  the  accuracy  of  state  estimation  and  the  speed  at  which  this 

accuracy  is  attained.  Examining  die  azimuth  error  estimation  curves  in 

Appendix  B,  we  see  that  in  every  case  the  incremental  velocity  (DVKAL) 

estimator  is  si  perior  or  equal  in  both  accuracy  attained  and  speed  of 

response,  Traisients  that  appear  in  the  individual  runs  are  due  to  system 

dynamics.  The  estimates  and  tn'i  values  for  the  azimuth  error  at  the  end 

of  the  alignment  period  (300  seconds)  as  well  as  the  R1C  estimate  errors 

are  given  in  Table  3.  If  we  define  the  acceptable  "criterion-of -goodness" 

as  estimating  the  azimuth  error  to  within  five  minutes  of  arc,  we  see 

that  adequate  performance  is  attained  on  only  three  of  the  TIKAL  runs. 

This  also  is  illustrated  by  the  RHS  error  term  which  is  greater  than  the 

required  bound.  In  comparing  the  response  tises,  we  observe  from 

Figures  12  and  13  that  the  average  square  of  the  DVKAL  azimuth  estimation 

2 

error  is  within  the  desired  bound  (25  min  )  after  approximately  200  seconds 


TIKHL  LEVEL  TILT 
XU)  ERROR  3QURREO 
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Figure  11.  DVKAL  level  tilt  covoriance  and  error  squared 
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Figure  12.  T1KAI.  azimuth  f'rror  covariance  and  error  squared 
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Table  3.  Azimuth  error  (oiautes  of  arc)  after  300  seconds 


Run 

True  Value 

DVKAL  Estimate 

TIKAL  Estimate 

1 

47.208 

45.239 

57.320 

2 

89.544 

87.214 

91.554 

3 

-5.046 

-5.262 

-8.219 

4 

33.24C 

32.650 

39.217 

5 

-47.490 

-47.280 

-54.528 

6 

-81.720 

-82.121 

-83.994 

7 

13.332 

12.240 

4.174 

8 

42.780 

40.228 

33.007 

9 

58.116 

54.924 

64.869 

10 

-106.248 

-104.198 

-112.373 

RMS  Error 

1.789 

6.856 

while  the  TIKAL  error  never  reaches  an  acceptable  level. 

Results  for  the  level  tilts  aie  similar  except  that  the  estimate 
is  obtained  more  rapidly.  The  reason  for  this  is  that  the  only  measure 
that  the  Kalman  filter  has  of  aziwch  error  is  its  connection  to  the  slope 
of  the  level  tilt.  Thus  we  expect  that  the  filter  must  first  estimate 
the  tilt  correctly  before  it  is  able  to  accurately  estimate  the  azimuth 
error  (which  is  constant). 

The  above  results  are  not  surprising  as  the  TIKAL  measurement  algorithm, 
in  effect,  models  the  physical  situation  incorrectly.  The  signal  physically 
available  to  the  TIKAL  estimator  is  a  measure  of  the  total  integrated 
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acceleration  from  the  accelerometer  (Figure  8)  and  is  dependent  upon 
x4(t)  and  x^(0).  Since  the  estimator  has  no  a  priori  knowledge  of  x^, 
its  best  estimate  of  x^(0)  is  zero.  Thus  we  suspect  that  as  x^(0) 
becomes  larger  the  performance  suffers  accordingly.  The  DVKAL  scheme 
is  not  affected  by  this  problem  as  it  uses  an  incremental  measurement 
and  estimates  the  necessary  present  and  delayed  states. 

Another  source  of  error  is  the  "residual  initial  condition"  which 
is  the  random  initial  voltage  on  the  accelerometer  integrator.  We  see 
from  Figure  9  that  the  resulting  integrated  acceleration  function  is  not 
zero  initially  as  is  assumed  in  the  TIKAL  estimator.  Since  the  DVKAL 
technique  ignores  the  first  AV  measurement,  errors  due  to  this  residual 
quantity  are  also  avoided. 

The  random  values  used  for  x^(0)  and  the  ,lresidual  initial  condition" 
as  well  as  the  final  squared  difference  between  the  azimuth  error  true 
and  estimated  values  are  listed  in  Table  4.  The  actual  numbers  listed 
for  x^(0)  and  the  "residual  initial  condition"  are  unity  variance  random 
numbers  with  the  necessary  additional  scaling  shown  at  the  top  of  the  table 
(the  "residual  initial  condition"  does  not  have  a  truly  normal  distribu¬ 
tion  as  it  is  not  allowed  to  exceed  tsie  granularity  AV).  As  seen  from 
Table  A  there  is  a  strong  correlation  between  large  final  errors  in  the 
TIKAL  estimation  and  high  initial  values  of  x^,  thus  verifying  our 
previous  suspicions.  This  effect  is  illustrated  by  taking  simulation 
Run  1  (which  had  a  large  value  for  x^(0)  )  and  setting  x^(0)  =  0.  The 
result  is  shown  in  Figure  14.  By  comparing  Figure  14  with  the  corre¬ 
sponding  graph  in  Appendix  B,  the  dramatic  improvement  in  performance 
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Table  4.  Final  Xj  simulation  errors  (ainctes  )  and  selected 
initial  conditions 


Run 

DVKAL 

Error  Squared 

TIKAL 

Error  Squared 

vo) 

(X0.69) 

Residual  Initial 

Condition  (x  AV) 

1 

3.876 

102.243 

-1.215 

-0.384 

2 

5.429 

4.039 

-0.016 

0.490 

3 

0.047 

10.067 

-0.153 

-0.735 

4 

0.348 

35.730 

-0.644 

0.088 

5 

0.044 

49.533 

0.827 

-0.064 

6 

0.161 

5.173 

0.343 

0.042 

7 

1.192 

83.860 

0.894 

-0.195 

8 

6.511 

95.516 

1.317 

0.707 

9 

10.137 

45.603 

-0.191 

0.969 

10 

4.202 

37.513 

0.312 

-0.682 

is  evident.  Note  that  the  DVKAL  performance  is  essentially  unchanged 
thus  supporting  our  earlier  allegation  that  the  DVKAL  estimator  is 
primarily  independent  of  x^(0).  Similar  results  are  observed  by  setting 
the  "residual  initial  condition”  equal  to  zero  in  runs  which  originally 
had  large  values  for  that  quantity. 

A  possible  source  of  error  in  the  DVKAL  estimator  arises  due  to  the 
measurement  noise  which  is  added  to  the  routine  to  insure  a  "safe” 

measurement  model,  it  is  possible  for  this  added  noise  to  slow  the  rate  at 
which  the  covariance  terms  decrease  and  thus  increase  the  uncertainty 
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associated  with  the  state  estimates.  In  (15./  and  (10),  the  added  noise 
decreases  the  gain,  b^,  thereby  weighting  the  measurement  term  less. 
Improvement  in  DVKAL  performance  for  Run  9  is  illustrated  in  Figure  15 
by  setting  the  measurement  noise  equal  to  zero.  The  D7KAL  response  with 
noise  equal  to  10%  of  AV  is  included  for  comparison.  Thus,  at  least 
for  one  run,  improved  performance  does  result  from  decreased  measurement 
noise.  This  means  that  a  smaller  value  of  measurement  noise  should  have 
been  chosen  for  the  simulation.  A  study  was  not  conducted  to  determine 
optimum  magnitudes  for  this  noise  and  the  value  selected  was  done  so  on 
a  purely  heuristic  basis. 

From  Figures  12  and  13  it  is  observed  that  the  TIKAL  covariance 
quantities  approach  zero  faster  than  the  corresponding  DVXAL  values. 

Tnis  is  primarily  due  to  the  aforementioned  DVKAL  measurement  noise  and 
(for  the  process  selected,  at  least)  a  larger  number  of  TIKAL  measurements. 
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VI.  CONCLUSIONS 

Under  the  conditions  end  limitations  assumed  {e.g. ,  zero  accelerom¬ 
eter  biases),  the  incremental  velocity  (DVKAL)  measurement  algorithm 
is  superior  to  the  time  interval  (TIXALJ  sampling  scheme  in  both  accuracy 
of  estimation  and  speed  of  response.  Although  the  results  presented  have 
been  obtained  using  a  single  physical  process,  the  positive  nature  of 
these  results  indicates  the  need  for  a  more  extensive  study  to  further 
establish  the  performance  level  and  limitations  associated  with  the 
incremental  velocity  technique. 

The  errors  incurred  by  the  time  interval  algorithm  are  derived 
primarily  from  its  improper  modeling  of  the  physical  situation. 


VII.  LITERATURE  CITED 


t,L\  «,  *  »•«  **.‘***J'W . 


rv— •  V ‘ 


50 

VIII.  ACKNOWLEDGMENTS 

The  author  is  indebted  to  Dr.  R.  G.  Brown  for  suggesting  the  topic 
for  this  study  and  for  the  guidance  given  during  its  investigation. 

This  work  was  supported  by  the  Engineering  Research  Institute  through 
funds  provided  by  the  Office  of  Naval  Research,  contract  number  N00014- 
68A-0162. 


51 


IX.  APPENDIX  A:  SIMULATION  PROGRAM 


A.  Main  Program 
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APPENDIX  B:  GRAPHICAL  RESULTS 


A.  Azimuth  Error 
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Figure  21.  Run  6  -  azimuth  error 
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Figure  23 
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Figure  24,  Hun  9  -  azimuth  error 
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.  Level  Tilt 
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Figure  26.  Run  1  -  level  tilt 
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”  . “  Figure  27.  Run  2  -  level  tilt 


Figure  28.  K'in  1  -  level  tilt 


Figure  29.  Run  4  -  level  tilt 
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Figure  31 


Run  7  -  level 
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Figure  34.  Run  9  -  level  tilt 
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